1 About

License:

CC-BY

To cite:

Anderson, B, (2024) New Zealand generation share trends: Electricity Authority data. https://dataknut.github.io/gridCarbon/nzGenMixTrends.html

Code:

https://github.com/dataknut/gridCarbon

2 Introduction

Inspired by https://twitter.com/DrSimEvans/status/1508409309775994892

3 Data

This analysis uses the New Zealand Electricity Authority’s Generation data.

Values are in kWh per half hour https://www.emi.ea.govt.nz/Wholesale/Datasets/Generation/Generation_MD

As far as we can work out this data does not include distributed (i.e. non-grid connected) generation such as small scale wind, solar, hydro, biomass etc which is connected to the LV network. This means the EA data is likely to underestimate total generation and potentially underestimate the proportion of total generation that is renewable. It is possible that this could be fixed using embedded wind & solar generation data from the metered embedded generation data although this may still not include household level (micro)generation.

rmdParams$olderThan <- 7

If the data we have previously downloaded is more than 7` days old, re-download.

rmdParams$dataPath <- "~/Dropbox/data/NZ_ElecAuth/processed/"

# load all the files we have using rbindlist
# probably ought to do some filtering here
filesToLoad <- list.files(rmdParams$dataPath, pattern = ".csv", 
                          full.names = TRUE)

orig_DT <- rbindlist(lapply(filesToLoad, fread))
setnames(orig_DT, "rDateTime", "dv_dateTime")

message("Original data range from: ", min(orig_DT$dv_dateTime, na.rm = TRUE))

## Original data range from: 2020-01-01 00:15:00

message("...to: ", max(orig_DT$dv_dateTime , na.rm = TRUE))

## ...to: 2024-07-31 23:45:00

# make wide for ease of aggregation ----
# sum over fuel code
nzGenMix_wide_dt <- dcast(orig_DT[!is.na(dv_dateTime)], dv_dateTime ~ Fuel_Code, value.var = "kWh", fun.aggregate = sum)

# add derived variables used later ----
addDerivedVars <- function(dt){
  dt[, dv_year := lubridate::year(dv_dateTime)]
dt[, dv_date := lubridate::as_date(dv_dateTime)]
dt[, dv_month := lubridate::month(dv_dateTime)]
dt[, dv_hour := lubridate::hour(dv_dateTime)]
dt[, dv_hms := hms::as_hms(dv_dateTime)]
dt <- gridCarbon::add_season(dt, 
                                 dateVar = "dv_dateTime",
                                 h = "S")
dt <- gridCarbon::add_peakPeriod(dt, 
                                      dateTime = "dv_dateTime")
  return(dt)
}

nzGenMix_wide_dt <- addDerivedVars(nzGenMix_wide_dt)
nzGenMix_dt <- addDerivedVars(orig_DT[!is.na(dv_dateTime)])

#message("Remove incomplete years to avoid weird things in plots.")

message("Filtered data range from: ", min(nzGenMix_wide_dt$dv_dateTime))

## Filtered data range from: 2020-01-01 00:15:00

message("...to: ", max(nzGenMix_wide_dt$dv_dateTime))

## ...to: 2024-07-31 23:45:00

rmdParams$plotCap <- paste0("Data: NZ EA Generation ",
                   min(nzGenMix_wide_dt$dv_date), " - ",
                   max(nzGenMix_wide_dt$dv_date),
                            "\nPlot: @dataknut \nCode: https://github.com/dataknut/gridCarbon")

# > recent date cut ----
recentDateCut <- max(nzGenMix_wide_dt$dv_date - (6*30)) # roughly 6 months

Note:

the data may cover more years than we need

the data may contain partial years - BEWARE incomplete years in plots using annual totals or means across all months.

4 Recreating DrSimEvans’ plot

Inspired by https://twitter.com/DrSimEvans/status/1508409309775994892

Which explored the relationship between % renewable and fossil-fuel generation for Great Britain (see also https://dataknut.github.io/gridCarbon/gbGenMixTrends.html)

This looks like daily data.

Solar + wind (% of gen) vs coal + gas (and % of gen))

# add together the % and totals we want (half-hourly)

nzGenMix_wide_dt[, dv_total := Coal + Diesel + Gas + Geo + Hydro + Solar + Wind + Wood]
nzGenMix_wide_dt[, dv_coal_gas := Coal + Gas]
nzGenMix_wide_dt[, dv_coal_gas_pc := 100 * dv_coal_gas/dv_total]
nzGenMix_wide_dt[, dv_solar_wind := Solar + Wind]
nzGenMix_wide_dt[, dv_solar_wind_pc := 100 * dv_solar_wind/dv_total]

# keep the vars we want for clarity
temp <- nzGenMix_wide_dt[, .(dv_dateTime, dv_coal_gas, dv_solar_wind,
                    dv_coal_gas_pc, dv_solar_wind_pc, dv_total)]

temp[, dv_date := lubridate::date(dv_dateTime)]

# aggregate to daily data for plotting
plotDT <- temp[,
               .(mean_dv_solar_wind_pc = mean(dv_solar_wind_pc),
                 mean_dv_coal_gas_pc = mean(dv_coal_gas_pc),
                 total_dv_coal_gas = sum(dv_coal_gas),
                 total_dv_solar_wind = sum(dv_solar_wind),
                 total_dv_total = sum(dv_total),
                 nObs = .N), # to check for days with < 48 half hours
               keyby = .(dv_date)
               ]
plotDT[, dv_year := lubridate::year(dv_date)] # for plots
plotDT[, total_dv_coal_gas_pc := 100 * total_dv_coal_gas/total_dv_total] # daily %
plotDT[, total_dv_solar_wind_pc := 100 * total_dv_solar_wind/total_dv_total]

message("Check for days with less than 48 hours - this will be truncated data due to DST breaks. We hate DST breaks")

## Check for days with less than 48 hours - this will be truncated data due to DST breaks. We hate DST breaks

table(plotDT$nObs)

## 
##   46   48 
##    4 1670

Figure 4.1 shows the mean half-hourly % generation by each type per day. This is slightly convoluted - it is the mean of the sum of the 48 daily half-hourly values. Unfold the code above for clarity.

The smoothed curves are estimated for each year. The lines terminate at the maximum value for the year. I’m still trying to decide if they tell us anything useful.

ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = mean_dv_solar_wind_pc, 
                            y = mean_dv_coal_gas_pc,
                            colour = as.factor(dv_year),
                            alpha = dv_year)) +
  geom_point() +
  geom_smooth() +
  scale_colour_viridis_d(name = "Year") +
  guides(alpha = "none") +
  labs(x = "Solar & wind (mean % of half-hourly generation per day)",
       y = "Coal & gas (mean % of half-hourly generation per day)",
       caption = rmdParams$plotCap)

## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 4.1: Mean half-hourly % generation by each type per day

# save it
ggplot2::ggsave(filename = "meanNZrenewablesVsfossilHalfHourPC.png", 
                path = rmdParams$plotPath,
                height = 5)

## Saving 7 x 5 in image
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 4.2 shows the percentage of daily generation by type. This is less convoluted as it is the sum of generation per day for the two categories (solar + wind vs gas + coal) as a % of total daily generation.

Again the smoothed curve is estimated for each year.

ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = total_dv_solar_wind_pc, 
                            y = total_dv_coal_gas_pc,
                            colour = as.factor(dv_year),
                            alpha = dv_year)) +
  geom_point() +
  geom_smooth() +
  scale_colour_viridis_d(name = "Year") +
  guides(alpha = "none") +
  labs(x = "Solar & wind (% of total daily generation)",
       y = "Coal & gas (% of total daily generation)",
       caption = rmdParams$plotCap)

## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 4.2: Percentage of daily generation by type

ggplot2::ggsave(filename = "dailyNZpcGenMix.png", 
                path = rmdParams$plotPath,
                height = 5)

## Saving 7 x 5 in image
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

4.1 Half-hourly versions of the plot

Just cos we can… helpfully split into ‘peak’ and ‘off peak’ periods.

Peak period definitions:

Morning 07:00 - 09:00
Daytime 09:00 - 16:00
Evening 16:00 - 21:00
Night - all other times

Again the smoothed curve is estimated for each year (and demand period).

ggplot2::ggplot(nzGenMix_wide_dt[dv_year > 2011], aes(x = dv_solar_wind_pc, 
                            y = dv_coal_gas_pc,
                            alpha = dv_year,
                            colour = as.factor(dv_year))) +
  geom_point() +
  facet_wrap(. ~ dv_peakPeriod) +
  geom_smooth() +
  scale_colour_viridis_d(name = "Year") +
  guides(alpha = "none") +
  labs(x = "Solar & wind (% of half-hourly generation)",
       y = "Coal & gas (% of half-hourly generation)",
       caption = rmdParams$plotCap)

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 4.3: Percentage of half-hourly generation by type

ggplot2::ggsave(filename = "halfHourlyNZ_PCgenByPeakPeriod.png", 
                path = rmdParams$plotPath,
                height = 5)

## Saving 7 x 5 in image
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

5 Daily generation trends

plotDT <- nzGenMix_dt[, 
                  .(dailyMWh = sum(as.numeric(kWh/1000))), # MWh -> MWh
                  keyby = .(dv_year, dv_month, dv_date, dv_season, Fuel_Code)]

# annual sum for cross-check 
t <- plotDT[, .(sum_TWh = sum(dailyMWh/1000000),
                nMonths = uniqueN(dv_month)), keyby = .(Year = dv_year)]
t[, Year := as.factor(Year)]

knitr::kable(t, caption = "Total TWh per year for cross-check - beware partial years")

Table 5.1: Total TWh per year for cross-check - beware partial years
Year	sum_TWh	nMonths
2020	40.61776	12
2021	40.97438	12
2022	40.66525	12
2023	40.59371	12
2024	24.34475	7

# double check

Figure 5.1 shows individual trends for daily generation. Note that the y scale is adapted to each fuel source for clarity of trends.

p <- ggplot2::ggplot(plotDT, aes(x = dv_date, y = dailyMWh/1000, colour = Fuel_Code)) +
  geom_line() +
  facet_grid(Fuel_Code ~ ., scales="free_y") +
  scale_color_viridis_d(name = "Source") +
  labs(x = "Date",
       y = "Daily GWh",
       caption = rmdParams$plotCap)

p <- fix_date_axis(p)
p

Figure 5.1: Trends in daily GWh generation

ggplot2::ggsave(p,
                filename = "dailyNZGenBySource.png", 
                path = rmdParams$plotPath,
                height = 5)

## Saving 7 x 5 in image

Figure 5.2 is interactive (hover mouse).

plotly::ggplotly(p)

Figure 5.2: Trends in daily GWh generation

Figure 5.3 stacks them to give a sense of proportion…

p <- ggplot2::ggplot(plotDT, aes(x = dv_date, y = dailyMWh/1000, fill = Fuel_Code)) +
  geom_col(position = "stack") +
  scale_fill_viridis_d(name = "Source") +
  labs(x = "Date",
       y = "Daily GWh",
       caption = rmdParams$plotCap)

p <- fix_date_axis(p)
p

Figure 5.3: Trends in daily GWh generation (stacked)

ggplot2::ggsave(p,
                filename = "genNZTrendStackBySource.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

Figure 5.4 shows the % daily contribution

plotDT[, pcOfDaily := 100 * dailyMWh/sum(dailyMWh), keyby = .(dv_date)]

p <- ggplot2::ggplot(plotDT, aes(x = dv_date, y = round(pcOfDaily,2), colour = Fuel_Code)) +
  geom_line() +
  scale_colour_viridis_d(name = "Source") +
  labs(x = "Date",
       y = "% of total daily GWh",
       caption = rmdParams$plotCap) 

p <- fix_date_axis(p)
p

Figure 5.4: Trends in % daily contribution

ggplot2::ggsave(p,
                filename = "genNZTrendPCBySource.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

Figure 5.5 - interactive version

plotly::ggplotly(p)

Figure 5.5: Interactive version

Figure 5.6 - interactive version

rp <- ggplot2::ggplot(plotDT[dv_date > recentDateCut], aes(x = dv_date, y = round(pcOfDaily,2), colour = Fuel_Code)) +
  geom_line() +
  scale_colour_viridis_d(name = "Source") +
  labs(x = "Date",
       y = "% of total daily GWh",
       caption = rmdParams$plotCap) 

rp <- fix_date_axis(rp)
rp <- rp +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

plotly::ggplotly(rp)

Figure 5.6: Interactive version

6 Detailed trends

6.1 Renewable %

Defined how?

Renewable = wind + solar + hydro + geothermal + wood - but note that geothermal is not necessarily low carbon.

nzGenMix_wide_dt[, RENEWABLE := Geo + Hydro + Solar + Wind + Wood]

nzGenMix_wide_dt[, RENEWABLE_perc := 100*RENEWABLE/dv_total]

make_NgesoPlots <- function(dt, 
                            var = "RENEWABLE_perc", # defaults
                            lowColour = "green", 
                            highColour = "red", 
                            scaleLab = "Change me!", 
                            unit_lab = "??",
                            minMax = "max"){
  res <- list() # for the results
  # tile plot ----
  res$tile <- ggplot2::ggplot(dt, aes(x = dv_date, 
                                      y = dv_hms, fill = get(var))) +
    geom_tile() +
    theme(legend.position = "bottom") +
    scale_x_date(date_labels="%b %Y",date_breaks  ="3 month") +
    theme(axis.text.x = element_text(angle = 60, hjust = 1)) +
    scale_fill_continuous(name = paste0(scaleLab, " (", unit_lab, ")"),
                          high = highColour,
                          low = lowColour) +
    labs(x = "Date",
         y = "Half-hour",
         caption = rmdParams$plotCap)
  
  # if(minMax = "max"){
  #   res$tile <- res$tile + 
  #     scale_fill_continuous(name = paste0("Maximum ", scaleLab)
  #                           )
  # }
  #   if(minMax = "min"){
  #   res$tile <- res$tile + 
  #     scale_fill_continuous(name = paste0("Half-hourly ", scaleLab)
  #                           )
  #   }
  
  # line plots ----
  plotDT <- dt[,.(max = max(get(var)),
                  mean = mean(get(var)),
                  min = min(get(var))
  ), keyby = .(dv_date)]
  
  #> mean half-hourly ----
  yLab <- paste0("Mean half-hourly ", scaleLab, " (", unit_lab, ")")
  yh <- max(plotDT$mean)
  res$line_mean <- ggplot2::ggplot(plotDT, aes(x = dv_date, 
                                              y = mean, colour = mean)) +
    theme(legend.position = "bottom") +
    scale_x_date(date_labels="%b %Y",date_breaks  ="3 month") +
    theme(axis.text.x = element_text(angle = 60, hjust = 1)) +
    geom_line() +
    scale_color_continuous(name = yLab,
                             high = highColour,
                             low = lowColour) +
    labs(x = "Date",
         y = yLab,
         caption = rmdParams$plotCap) + 
      geom_smooth(aes(y = mean), colour = "grey") +
      geom_hline(yintercept = yh, colour = highColour) +
      annotate("text", x = mean(plotDT$dv_date), 
                      y = yh*1.05, 
               label = paste0("Maximum value: ", round(yh,1), unit_lab))
  
  #> max half-hourly ----
  yLab <- paste0("Maximum half-hourly ", scaleLab, " (", unit_lab, ")")
  yh <- max(plotDT$max)
  res$line_max <- ggplot2::ggplot(plotDT, aes(x = dv_date, y = max, colour = max)) +
    theme(legend.position = "bottom") +
    scale_x_date(date_labels="%b %Y",date_breaks  ="3 month") +
    theme(axis.text.x = element_text(angle = 60, hjust = 1)) +
    geom_line() +
    scale_color_continuous(name = yLab,
                             high = highColour,
                             low = lowColour) +
    labs(x = "Date",
         y = yLab,
         caption = rmdParams$plotCap) + 
      geom_smooth(aes(y = max), colour = "grey") +
      geom_hline(yintercept = yh, colour = highColour) +
      annotate("text", x = mean(plotDT$dv_date), 
                      y = yh*1.05, 
               label = paste0("Historical maximum: ", round(yh,1), unit_lab))
     
  #> min half-hourly ----
  yLab <- paste0("Minimum half-hourly ", scaleLab, " (", unit_lab, ")")
  yh <- min(plotDT$min)
  res$line_min <- ggplot2::ggplot(plotDT, aes(x = dv_date, y = min, colour = min)) +
    theme(legend.position = "bottom") +
    scale_x_date(date_labels="%b %Y",date_breaks  ="3 month") +
    theme(axis.text.x = element_text(angle = 60, hjust = 1)) +
    geom_line() +
    scale_color_continuous(name = yLab,
                             high = highColour,
                             low = lowColour) +
    labs(x = "Date",
         y = yLab,
         caption = rmdParams$plotCap) + 
      geom_smooth(aes(y = min), colour = "grey") +
      geom_hline(yintercept = yh, colour = highColour) +
      annotate("text", x = mean(plotDT$dv_date), 
                      y = yh*1.05, 
               label = paste0("Historical minimum: ", round(yh,1), unit_lab)
               )
  
  minDateTime <- dt[get(var) == min(plotDT$min), dv_dateTime]
  res$minDateReport <- paste0("Minimum: ", 
                              round(plotDT[min == min(plotDT$min), min],1), unit_lab,
                                    " at ", minDateTime)
  maxDateTime <- dt[get(var) == max(plotDT$max), dv_dateTime]
  res$maxDateReport <- paste0("Maximum: ", 
                              round(plotDT[max == max(plotDT$max), max],1), unit_lab,
                                    " at ", maxDateTime)

  
  return(res)
}

6.1.1 Renewable % (all years)

results <-  make_NgesoPlots(nzGenMix_wide_dt, 
                       var = "RENEWABLE_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Renewable generation",
                       unit = "%",
                       minMax = "max" # max line or min line?
)
results$tile

Figure 6.1: Half hourly % renewables, all years

results$line_mean

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.2: Mean half-hourly % renewables per day, all years

results$line_max

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.3: Maximum half-hourly % renewables per day, all years

results$maxDateReport

## [1] "Maximum: 99% at 2023-05-24 00:45:00"

results$line_min

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.4: Minimum % renewables, all years

results$minDateReport

## [1] "Minimum: 64.3% at 2021-05-02 23:15:00"

6.1.2 Renewable % (last 6 months)

results <-  make_NgesoPlots(nzGenMix_wide_dt[dv_date > recentDateCut], 
                       var = "RENEWABLE_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Renewable generation",
                       minMax = "max",
                       unit_lab = "%"
)
results$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.5: % renewables, last 6 months

results$line_mean +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.6: Mean % renewables, last 6 months

results$line_max+
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.7: Maximum % renewables, last 6 months

results$maxDateReport

## [1] "Maximum: 94.6% at 2024-02-05 10:45:00"

results$line_min +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.8: Minimum % renewables, last 6 months

results$minDateReport

## [1] "Minimum: 67.6% at 2024-07-23 01:15:00"

6.2 Wind %

6.2.1 Wind % (all years)

nzGenMix_wide_dt[, WIND_perc := 100*Wind/dv_total]
results <-  make_NgesoPlots(nzGenMix_wide_dt, 
                       var = "WIND_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Wind generation",
                       unit = "%",
                       minMax = "max" # max line or min line?
)
results$tile

Figure 6.9: % wind, all years

results$line_mean

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.10: Mean % wind per day, all years

results$line_max

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.11: Maximum % wind, all years

results$maxDateReport

## [1] "Maximum: 24.6% at 2024-04-26 04:15:00"

results$line_min

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.12: Minimum % wind, all years

results$minDateReport

## [1] "Minimum: 0% at 2020-10-29 10:45:00" "Minimum: 0% at 2020-10-29 11:15:00" "Minimum: 0% at 2020-10-29 11:45:00"
## [4] "Minimum: 0% at 2020-10-29 12:15:00" "Minimum: 0% at 2021-03-24 01:45:00" "Minimum: 0% at 2022-04-03 06:15:00"
## [7] "Minimum: 0% at 2022-04-03 06:45:00"

6.2.2 Wind % (last 6 months)

results <-  make_NgesoPlots(nzGenMix_wide_dt[dv_date > recentDateCut], 
                       var = "WIND_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Wind generation",
                       minMax = "max",
                       unit_lab = "%"
)
results$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.13: % wind, last 6 months

Figure 6.14 shows the mean half-hourly % of wind generation per day.

results$line_mean +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.14: Mean % wind, last 6 months

results$maxDateReport

## [1] "Maximum: 24.6% at 2024-04-26 04:15:00"

results$line_max+
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.15: Maximum % wind, last 6 months

results$maxDateReport

## [1] "Maximum: 24.6% at 2024-04-26 04:15:00"

results$line_min +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.16: Minimum % wind, last 6 months

results$minDateReport

## [1] "Minimum: 0% at 2024-04-01 12:45:00"

6.3 Solar %

6.3.1 Solar % (all years)

nzGenMix_wide_dt[, SOLAR_perc := 100*Solar/dv_total]
results <-  make_NgesoPlots(nzGenMix_wide_dt, 
                       var = "SOLAR_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Solar generation",
                       unit = "%",
                       minMax = "max" # max line or min line?
)
results$tile

Figure 6.17: % solar, all years

results$line_mean

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.18: Mean % solar, all years

results$line_max

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.19: Maximum % solar, all years

results$maxDateReport

## [1] "Maximum: 0.7% at 2024-07-28 15:15:00"

results$line_min

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

## Warning: Failed to fit group -1.
## Caused by error in `gam.reparam()`:
## ! NA/NaN/Inf in foreign function call (arg 3)

Figure 6.20: Minimum % solar, all years

#results$minDateReport

6.3.2 Solar % (last 6 months)

results <-  make_NgesoPlots(nzGenMix_wide_dt[dv_date > recentDateCut], 
                       var = "SOLAR_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Solar generation",
                       minMax = "max",
                       unit_lab = "%"
)
results$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.21: % solar, last 6 months

Figure 6.22 shows the mean daily % of generation which is (grid) solar.

results$line_mean+
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.22: Mean % solar, last 6 months

results$line_max+
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.23: Maximum % solar, last 6 months

results$maxDateReport

## [1] "Maximum: 0.7% at 2024-07-28 15:15:00"

results$line_min +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.24: Minimum % solar, last 6 months

#results$minDateReport

6.4 Hydro %

6.4.1 Hydro % (all years)

nzGenMix_wide_dt[, HYDRO_perc := 100*Hydro/dv_total]
results <-  make_NgesoPlots(nzGenMix_wide_dt, 
                       var = "HYDRO_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Hydro generation",
                       unit = "%",
                       minMax = "max" # max line or min line?
)
results$tile

Figure 6.25: % solar, all years

results$line_mean

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.26: Mean % hydro, all years

results$line_max

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.27: Maximum % hydro, all years

results$maxDateReport

## [1] "Maximum: 83.2% at 2022-12-03 11:45:00"

results$line_min

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.28: Minimum % hydro, all years

results$minDateReport

## [1] "Minimum: 25.8% at 2024-07-29 03:15:00"

6.4.2 Hydro % (last 6 months)

results <-  make_NgesoPlots(nzGenMix_wide_dt[dv_date > recentDateCut], 
                       var = "HYDRO_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Hydro generation",
                       minMax = "max",
                       unit_lab = "%"
)
results$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.29: % solar, last 6 months

results$line_mean +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.30: Mean % hydro, last 6 months

results$line_max+
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.31: Maximum % hydro, last 6 months

results$maxDateReport

## [1] "Maximum: 70.8% at 2024-02-07 17:15:00"

results$line_min +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.32: Minimum % hydro, last 6 months

results$minDateReport

## [1] "Minimum: 25.8% at 2024-07-29 03:15:00"

6.5 Geothermal %

6.5.1 Geothermal % (all years)

nzGenMix_wide_dt[, GEO_perc := 100*Geo/dv_total]
results <-  make_NgesoPlots(nzGenMix_wide_dt, 
                       var = "GEO_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Geothermal generation",
                       unit = "%",
                       minMax = "max" # max line or min line?
)
results$tile

Figure 6.33: % solar, all years

results$line_mean

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.34: Mean % geothermal, all years

results$line_max

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.35: Maximum % geothermal, all years

results$maxDateReport

## [1] "Maximum: 34.2% at 2020-04-13 03:45:00"

results$line_min

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.36: Minimum % geothermal, all years

results$minDateReport

## [1] "Minimum: 9.6% at 2023-07-07 17:45:00"

6.5.2 Geothermal % (last 6 months)

results <-  make_NgesoPlots(nzGenMix_wide_dt[dv_date > recentDateCut], 
                       var = "GEO_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Geothermal generation",
                       minMax = "max",
                       unit_lab = "%"
)
results$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.37: % geothermal, last 6 months

results$line_mean +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.38: Mean % geothermal, last 6 months

results$line_max+
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.39: Maximum % geothermal, last 6 months

results$maxDateReport

## [1] "Maximum: 31.4% at 2024-06-02 03:45:00"

results$line_min +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.40: Minimum % geothermal, last 6 months

results$minDateReport

## [1] "Minimum: 12.5% at 2024-06-24 18:15:00"

6.6 Gas %

6.6.1 Gas % (all years)

nzGenMix_wide_dt[, GAS_perc := 100*Gas/dv_total]
results <-  make_NgesoPlots(nzGenMix_wide_dt, 
                       var = "GAS_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Gas generation",
                       unit = "%",
                       minMax = "max" # max line or min line?
)
results$tile

Figure 6.41: % gas, all years

results$line_mean

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.42: Mean % gas, all years

results$line_max

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.43: Maximum % gas, all years

results$maxDateReport

## [1] "Maximum: 25.1% at 2020-07-14 23:45:00"

results$line_min

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.44: Minimum % gas, all years

results$minDateReport

## [1] "Minimum: 0.9% at 2023-11-06 14:45:00"

6.6.2 Gas % (last 6 months)

results <-  make_NgesoPlots(nzGenMix_wide_dt[dv_date > recentDateCut], 
                       var = "GAS_perc",
                       lowColour = "green",
                       highColour = "red",
                       scaleLab = "Gas generation",
                       minMax = "max",
                       unit_lab = "%"
)
results$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.45: % geothermal, last 6 months

results$line_mean +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.46: Mean % gas, last 6 months

results$line_max+
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.47: Maximum % gas, last 6 months

results$maxDateReport

## [1] "Maximum: 19.9% at 2024-06-23 03:15:00"

results$line_min +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.48: Minimum % gas, last 6 months

results$minDateReport

## [1] "Minimum: 1.9% at 2024-04-11 00:45:00"

6.7 Coal %

6.7.1 Coal % (all years)

nzGenMix_wide_dt[, COAL_perc := 100*Coal/dv_total]
results <-  make_NgesoPlots(nzGenMix_wide_dt, 
                       var = "COAL_perc",
                       lowColour = "green",
                       highColour = "red",
                       scaleLab = "Coal generation",
                       unit = "%",
                       minMax = "max" # max line or min line?
)
results$tile+
  scale_x_date(date_labels="%b-%Y")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.49: % coal, all years

results$line_mean +
  scale_x_date(date_labels="%b-%Y")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.50: Mean % coal, all years

results$line_max+
  scale_x_date(date_labels="%b-%Y")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.51: Maximum % coal, all years

results$maxDateReport

## [1] "Maximum: 20.6% at 2021-04-21 03:45:00"

results$line_min +
  scale_x_date(date_labels="%b-%Y")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.52: Minimum % coal, all years

#results$minDateReport

6.7.2 Coal % (last 6 months)

results <-  make_NgesoPlots(nzGenMix_wide_dt[dv_date > recentDateCut], 
                       var = "COAL_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Coal generation",
                       minMax = "max",
                       unit_lab = "%"
)
results$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.53: % coal, last 6 months

results$line_mean +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.54: Mean % coal, last 6 months

results$line_max+
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.55: Maximum % coal, last 6 months

results$maxDateReport

## [1] "Maximum: 16.4% at 2024-07-26 01:15:00"

results$line_min +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="1 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 6.56: Minimum % coal, last 6 months

#results$minDateReport

6.8 All ‘fuels’ half-hourly trends

These plots are a denser version of (e.g.) Figure 5.1 because they plot half-hourly %.

# calculate % contribution in each half hour
# long form data
plotDT <- nzGenMix_dt[, .(GWh = sum(kWh/1000000, na.rm = TRUE)), 
                      keyby = .(Fuel_Code, dv_dateTime)]
plotDT[, pc := 100*GWh/sum(GWh), keyby = .(dv_dateTime)]

p <- ggplot2::ggplot(plotDT,
                aes(x = dv_dateTime,
                    y = pc,
                    colour = Fuel_Code)) +
  geom_line() +
  scale_colour_viridis_d() +
  labs(y = "%",
       x = "Half hour")

fix_dateTime_axis(p)

Figure 6.57: Half-hourly % contribution by fuel over time (all years)

p <- ggplot2::ggplot(plotDT,
                aes(x = dv_dateTime,
                    y = pc,
                    fill = Fuel_Code)) +
  geom_col(position = "stack") +
  scale_fill_viridis_d() +
  labs(y = "%",
       x = "Half hour")

fix_dateTime_axis(p) +
  coord_flip()

# calculate % contribution in each half hour
p <- ggplot2::ggplot(plotDT[dv_dateTime > recentDateCut],
                aes(x = dv_dateTime,
                    y = pc,
                    colour = Fuel_Code)) +
  geom_line() +
    scale_colour_viridis_d() +
  labs(y = "%",
       x = "Half hour")

fix_dateTime_axis(p, breaks = "1 month")

Figure 6.58: Half-hourly % contribution by fuel over time (last 6 months)

p <- ggplot2::ggplot(plotDT[dv_dateTime > recentDateCut],
                aes(x = dv_dateTime,
                    y = pc,
                    fill = Fuel_Code)) +
  geom_col(position = "stack") +
  scale_fill_viridis_d() +
  labs(y = "%",
       x = "Half hour")

fix_dateTime_axis(p, breaks = "1 month") +
  coord_flip()

6.8.1 Hydro vs wind vs grid solar

p <- ggplot2::ggplot(plotDT[Fuel_Code == "Wind" | Fuel_Code == "Hydro" |
                              Fuel_Code == "Solar"],
                aes(x = dv_dateTime,
                    y = pc,
                    colour = Fuel_Code)) +
  geom_line() +
  scale_colour_manual(values = c("lightblue", 
                                 "orange", 
                                 "green") # what colour is the wind?!
                      ) +
  labs(y = "% of generation",
       x = "Half hour")
fix_dateTime_axis(p)

Figure 6.59: Comparison of hydro vs wind vs solar % half-hourly contribution over time (all years)

Figure 6.60 shows the recent half-hourly % of generation from hydro, wind and grid solar.

p <- ggplot2::ggplot(plotDT[dv_dateTime > recentDateCut & (Fuel_Code == "Wind" | Fuel_Code == "Hydro" | Fuel_Code == "Solar")],
                aes(x = dv_dateTime,
                    y = pc,
                    colour = Fuel_Code)) +
  geom_line() +
  scale_colour_manual(values = c("lightblue", 
                                 "orange", 
                                 "green") # what colour is the wind?!
                      ) +
  labs(y = "% of generation",
       x = "Half hour")
fix_dateTime_axis(p, breaks = "1 month")

Figure 6.60: Comparison of hydro vs wind vs solar % half-hourly contribution over time (last 6 months)

nzGenMix_wide_dt[, ratio := WIND_perc/HYDRO_perc]
yMax <- max(nzGenMix_wide_dt$ratio)
p <- ggplot2::ggplot(nzGenMix_wide_dt,
                aes(x = dv_dateTime,
                    y =  100*ratio,
                    colour = dv_season,
                    group = dv_season)) +
  geom_point() +
  geom_hline(yintercept = 100*yMax) +
  annotate("text", x = mean(nzGenMix_wide_dt$dv_dateTime), 
                      y = 100 * yMax*1.05, 
               label = paste0("Maximum: ", 100*round(yMax,2), "%"))+
  labs(y = "Wind GWh as a % of Hydro GWh",
       x = "Half hour")
fix_dateTime_axis(p)

Figure 6.61: Wind GWh as a % of hydro GWh (all years)

yMax <- max(nzGenMix_wide_dt[dv_dateTime > recentDateCut]$ratio)
p <- ggplot2::ggplot(nzGenMix_wide_dt[dv_dateTime > recentDateCut],
                aes(x = dv_dateTime,
                    y =  100*ratio,
                    colour = dv_season,
                    group = dv_season)) +
  geom_point() +
  geom_hline(yintercept = 100*yMax) +
  annotate("text", x = mean(nzGenMix_wide_dt[dv_dateTime > recentDateCut]$dv_dateTime), 
                      y = 100 * yMax*1.05, 
               label = paste0("Maximum: ", 100*round(yMax,2), "%"))+
  labs(y = "Wind GWh as a % of Hydro GWh",
       x = "Half hour")
fix_dateTime_axis(p, breaks = "1 month")

Figure 6.62: Wind GWh as a % of hydro GWh (last 6 months)

7 Fuel types and time of day ‘demand’

Do we see a relationship between generation fuels and peak demand? This will be mediated by the way the electricity market works.

We may find wind curtailment (not visible here) at low demand periods where geothermal can’t be shut off.

First, what is the general average shape of overall generation?

‘Peak periods’ shown as shaded rectangles

# orginal values are kWh
plotDT <- nzGenMix_wide_dt[, .(mean_total_GWh = mean(dv_total/1000000),
                               mean_renewables_GWh = mean(RENEWABLE/1000000),
                               mean_renewables_pc = mean(RENEWABLE_perc),
                               mean_hydro_GWh = mean(Hydro/1000000),
                               mean_hydro_pc = mean(HYDRO_perc),
                               mean_coal_GWh = mean(Coal/1000000),
                               mean_coal_pc = mean(COAL_perc),
                               mean_gas_GWh = mean(Gas/1000000),
                               mean_gas_pc = mean(GAS_perc)), keyby = .(dv_hms, dv_year, dv_season, dv_peakPeriod)]

make_profilePlot <- function(dt, yVar = "mean_renewables_GWh", yVarLab = "change me"){
  plotDT <- dt[, .(x = dv_hms, y = get(yVar), alpha = dv_year, 
                   group = dv_year,
                   facet1 = dv_season,
                   colour = dv_peakPeriod)]
  p <- ggplot2::ggplot(plotDT, aes(x = x, 
                                   y = y,
                            alpha = alpha,
                            group = group)) +
    geom_rect(aes(xmin=hms::as_hms("07:00:00"), xmax=hms::as_hms("09:00:00"),
                  ymin=min(plotDT$y), ymax=max(plotDT$y)),
                              fill='pink', alpha=0.1) +
        geom_rect(aes(xmin=hms::as_hms("16:00:00"), xmax=hms::as_hms("21:00:00"),
                  ymin=min(plotDT$y), ymax=max(plotDT$y)),
                              fill='pink', alpha=0.1) +
  geom_line() +
  scale_alpha_continuous(name = "Year") +
  facet_wrap(facet1 ~ .) +
    labs(x = "Time of day",
       y = yVarLab,
       caption = rmdParams$plotCap) +
  theme(axis.text.x = element_text(angle = 45, hjust=1)) +
    scale_x_time(breaks = scales::breaks_width("120 min")) 
  return(p)
}

make_profilePlot(plotDT, "mean_total_GWh", "Mean total (GWh)")

Figure 7.1: Mean half-hourly generation by year and season

Contribution of renewables…

make_profilePlot(plotDT, "mean_renewables_GWh", "Mean renewables (GWh)")

Figure 7.2: Mean half-hourly % renewable generation by year and season

make_profilePlot(plotDT, "mean_renewables_pc", "Mean % renewables")

Figure 7.3: Mean half-hourly % renewable generation by year and season

What is hydro’s contribution?

make_profilePlot(plotDT, "mean_hydro_GWh", "Mean hydro (GWh)")

Figure 7.4: Mean half-hourly % hydro generation by year and season

make_profilePlot(plotDT, "mean_hydro_pc", "Mean % hydro")

Figure 7.5: Mean half-hourly % hydro generation by year and season

What is coal’s contribution?

make_profilePlot(plotDT, "mean_coal_GWh", "Mean coal (GWh)")

Figure 7.6: Mean half-hourly % coal generation by year and season

make_profilePlot(plotDT, "mean_coal_pc", "Mean % coal")

Figure 7.7: Mean half-hourly % coal generation by year and season

What is gas’s contribution?

make_profilePlot(plotDT, "mean_gas_GWh", "Mean gas (GWh")

Figure 7.8: Mean half-hourly % gas generation by year and season

make_profilePlot(plotDT, "mean_gas_pc", "Mean % gas")

Figure 7.9: Mean half-hourly % gas generation by year and season

Relationship between hydro and overall generation

ggplot2::ggplot(nzGenMix_wide_dt, aes(y = Hydro/1000000, x = dv_total/1000000, colour = dv_peakPeriod)) +
  geom_point() +
  facet_grid(dv_year ~ dv_season) +
  labs(y = "Hydro GWh per half hour",
       x = "Total GWh per half hour")

YMMV on 7.10

plotDT <- nzGenMix_wide_dt[, .(mean_renewables_GWh = mean(RENEWABLE/1000000),
                 mean_generation_GWh = mean(dv_total/1000000)),
             keyby = .(dv_year, dv_hms, dv_peakPeriod, dv_season)]

ggplot2::ggplot(plotDT, aes(x = mean_generation_GWh , y = mean_renewables_GWh,
                                         colour = dv_peakPeriod)) +
  geom_point() +
  scale_color_discrete(name = "Period") +
  facet_grid(dv_season ~ dv_year) +
  labs(x = "Mean half-hourly generation (GWh)",
       y = "Mean half-hourly renewables (GWh)",
       caption = rmdParams$plotCap)

Figure 7.10: Mean renewables vs Mean total generation

8 Summary

That’s it.

You might want to look at recent academic research on this topic:

(Staffell 2017)
(Staffell and Pfenninger 2018)

9 Annex

9.1 Carbon intensity trends

XX need to create a half-hourly CARBON_INTENSITY variable XX Requires: * per-fuel CI

Currently: https://environment.govt.nz/assets/publications/Measuring-emissions-guidance-August-2022/Emission-factors-workbook-Measuring-emissions-guidance-August-2022.xlsx

“The grid-average emission factor best reflects the carbon dioxide equivalent emissions associated with the generation of a unit of electricity purchased from the national grid in New Zealand in 2020”

“We calculate purchased electricity emission factors on a calendar-year basis and based on the average grid mix of generation types for calendar years. The emission factor accounts for the emissions from fuel combustion at thermal power stations and fugitive emissions from the generation of geothermal electricity. Thermal electricity is generated by burning fossil fuels.

The emission factor for purchased grid-average electricity does not include transmission and distribution losses. ”

“This emission factor also doesn’t reflect the real-world factors that influence the carbon intensity of the grid such as time of year, time of day and geographical area. Therefore, a grid -average emission factor may over or underestimate your organisation’s GHG emissions.”

2020 value = 0.120 kg CO2/kWh

Instead use a model derived from the GB MWh generation mix for 2023.

Coefficients:

Intercept: 41.231171
GAS_perc: 3.571570
COAL_perc: 10.050410
WIND_perc: -0.451285
HYDRO_perc: 1.147793 # why is hydro +ve?
SOLAR_perc: -0.335122

This is (obviously) completely broken as HYDRO_perc has a +ve coefficient, the model does not include a coefficient for Geothermal or Wood and the GB mix includes Nuclear, Interconnection and Storage. #YMMV

Need to fix :-)

But let’s imagine it is valid… what would we see?

# should be annual and fuel specific
# 2020 value

nzGenMix_wide_dt[, CARBON_INTENSITY := 41.231171 +
                   GAS_perc   *  3.571570   + 
                   COAL_perc   * 10.050410   +
                   WIND_perc  * -0.451285  +
                   HYDRO_perc *  1.147793 + # why is hydro +ve?
                   SOLAR_perc * -0.335122]

summary(nzGenMix_wide_dt$CARBON_INTENSITY)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   100.2   151.2   184.0   187.2   218.1   339.1

9.1.1 Overall trends

Seasons:

Summer: Jun - Aug
Autumn: Sept - Nov
Winter: Dec - Feb
Spring: Mar - May

plotDT <- nzGenMix_wide_dt[, .(mean_CI = mean(CARBON_INTENSITY)), keyby = .(dv_year, dv_date, dv_season)]

ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = dv_date, y = mean_CI, colour = dv_season)) +
  geom_point() +
  scale_color_viridis_d(name = "Season") +
  theme(legend.position = "bottom") +
  labs(x = "Date",
       y = "Mean g CO2e/kWh",
       caption = rmdParams$plotCap)

Figure 9.1: Mean half-hourly carbon intensity per day

ggplot2::ggsave(filename = "meanCiTrendBySeason.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

Re-draw as a boxplot of mean daily CI by month - plotted at month start. 9.2 shows outliers nicely.

plotDT <- nzGenMix_wide_dt[, .(mean_CI = mean(CARBON_INTENSITY)), 
                     keyby = .(dv_year, 
                               dv_date,
                               dv_season)]
plotDT[, dv_month := lubridate::floor_date(dv_date, unit = "months"),]
ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = dv_month, 
                                            y = mean_CI, 
                                            group = dv_month,
                                            colour = dv_season)) +
  geom_boxplot() +
  scale_color_viridis_d(name = "Season") +
  theme(legend.position = "bottom") +
  labs(x = "Month",
       y = "Mean g CO2e/kWh",
       caption = rmdParams$plotCap)

Figure 9.2: Mean daily CI boxplots per month

ggplot2::ggsave(filename = "meanCiTrendByMonthSeason.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

Figure 9.3 shows half-hourly carbon intensity over time.

#plotDT <- nzGenMix_wide_dt[, .()]
ci <- make_NgesoPlots(nzGenMix_wide_dt, 
                      var = "CARBON_INTENSITY",
                      lowColour = "green",
                      highColour = "brown",
                      scaleLab = "Carbon intensity (g CO2/MW)",
                      minMax = "min"
)
ci$tile

Figure 9.3: Half-hourly carbon intensity over time

ci$line

## NULL

ci <- make_NgesoPlots(nzGenMix_wide_dt[dv_year > 2021], 
                      var = "CARBON_INTENSITY",
                      lowColour = "green",
                      highColour = "brown",
                      scaleLab = "Carbon intensity (g CO2/MW)",
                      minMax = "min"
)
ci$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="3 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

ci$line +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="3 month")

## NULL

9.1.2 Annual by peak period

Mean carbon intensity per year and season within peak period.

Figure 9.4 suggests evening peak periods still have slightly higher carbon intensity and the shape of the reduction curves differ by season although rather less by period. Interestingly the sustained reduction in carbon intensity in Summer has leveled off.

plotDT <- nzGenMix_wide_dt[, .(mean_CI = mean(CARBON_INTENSITY)), keyby = .(dv_year, dv_peakPeriod, dv_season)]

ggplot2::ggplot(plotDT, aes(x = dv_year, y = mean_CI, colour = dv_peakPeriod)) +
  geom_line() +
  scale_color_viridis_d(name = "Peak period") +
  facet_wrap(. ~ dv_season) +
  labs(x = "Year",
       y = "Mean g CO2e/kWh",
       caption = rmdParams$plotCap)

Figure 9.4: Mean half-hourly carbon intensity by peak period and season 2012-2022

ggplot2::ggsave(filename = "annualMeanCIByPeak.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

9.2 Data descriptors

9.2.1 NZ gen mix data

skimr::skim(nzGenMix_wide_dt)

Table 9.1: Data summary
Name	nzGenMix_wide_dt
Number of rows	80344
Number of columns	32
Key	dv_dateTime
_______________________
Column type frequency:
Date	1
difftime	2
factor	2
numeric	26
POSIXct	1
________________________
Group variables	None

Variable type: Date

skim_variable	n_missing	complete_rate	min	max	median	n_unique
dv_date	0	1	2020-01-01	2024-07-31	2022-04-16	1674

Variable type: difftime

skim_variable	n_missing	complete_rate	min	max	median	n_unique
dv_hms	0	1	900 secs	85500 secs	42300 secs	48
hms	0	1	900 secs	85500 secs	42300 secs	48

Variable type: factor

skim_variable	n_missing	complete_rate	ordered	n_unique	top_counts
dv_season	0	1	FALSE	4	Aut: 22080, Win: 20592, Sum: 20208, Spr: 17464
dv_peakPeriod	0	1	FALSE	5	Ear: 23436, Day: 23436, Eve: 13392, Lat: 13384

Variable type: numeric

skim_variable	complete_rate	mean	sd	p0	p25	p50	p75	p100	hist
Coal	1	102551.40	95511.68	0.00	0.00	85433.00	174755.50	365201.00	▇▅▃▂▁
Diesel	1	350.50	3382.94	0.00	0.00	0.00	0.00	79213.00	▇▁▁▁▁
Gas	1	257958.57	124946.18	18133.00	171104.75	250195.00	339700.00	646409.00	▃▇▆▂▁
Geo	1	427833.33	35953.24	224851.00	410573.50	431223.50	447368.25	550997.00	▁▁▅▇▁
Hydro	1	1374752.44	346190.15	465067.00	1105091.00	1393511.50	1623827.75	2363701.00	▂▆▇▅▁
Solar	1	193.34	1345.51	0.00	0.00	0.00	0.00	19086.00	▇▁▁▁▁
Wind	1	153262.07	98054.21	0.00	70677.00	139026.00	223276.75	477396.00	▇▇▅▂▁
Wood	1	13027.73	5460.84	0.00	11880.00	15231.00	16670.00	19561.00	▂▂▁▆▇
dv_year	1	2021.82	1.34	2020.00	2021.00	2022.00	2023.00	2024.00	▇▇▇▇▅
dv_month	1	6.20	3.40	1.00	3.00	6.00	9.00	12.00	▇▆▆▅▆
dv_hour	1	11.50	6.92	0.00	5.00	11.00	17.00	23.00	▇▇▆▇▇
dv_total	1	2329929.38	404976.50	1324732.00	1990030.50	2361135.00	2608821.75	3626574.00	▂▆▇▂▁
dv_coal_gas	1	360509.97	183929.44	19460.00	212909.75	343146.50	495899.25	962020.00	▆▇▇▂▁
dv_coal_gas_pc	1	15.18	6.77	1.02	9.77	15.23	20.39	35.68	▃▇▇▅▁
dv_solar_wind	1	153455.41	98092.74	0.00	70841.75	139161.50	223397.00	477396.00	▇▇▅▂▁
dv_solar_wind_pc	1	6.77	4.50	0.00	3.08	5.98	9.79	24.57	▇▆▃▁▁
RENEWABLE	1	1969068.91	337707.31	1143537.00	1700021.50	1989244.00	2216382.00	3062302.00	▂▆▇▃▁
RENEWABLE_perc	1	84.81	6.79	64.32	79.60	84.77	90.23	98.98	▁▅▇▇▃
WIND_perc	1	6.77	4.50	0.00	3.07	5.97	9.79	24.57	▇▆▃▁▁
SOLAR_perc	1	0.01	0.05	0.00	0.00	0.00	0.00	0.74	▇▁▁▁▁
HYDRO_perc	1	58.50	8.22	25.83	53.36	58.74	63.94	83.17	▁▂▇▆▁
GEO_perc	1	18.95	3.84	9.56	16.06	18.14	21.67	34.20	▂▇▅▂▁
GAS_perc	1	10.91	4.63	0.88	8.00	11.01	14.10	25.10	▃▇▇▃▁
COAL_perc	1	4.27	3.93	0.00	0.00	3.64	7.03	20.65	▇▃▂▁▁
ratio	1	0.12	0.10	0.00	0.05	0.10	0.17	0.78	▇▃▁▁▁
CARBON_INTENSITY	1	187.18	42.25	100.15	151.21	184.01	218.09	339.08	▅▇▆▂▁

Variable type: POSIXct

skim_variable	n_missing	complete_rate	min	max	median	n_unique
dv_dateTime	0	1	2020-01-01 00:15:00	2024-07-31 23:45:00	2022-04-17	80344

9.3 R environment

Packages etc:

base R (R Core Team 2016)
bookdown (Xie 2016a)
data.table (Dowle et al. 2015)
ggplot2 (Wickham 2009)
here (Müller 2017)
hms (Müller 2018)
knitr (Xie 2016b)
lubridate (Grolemund and Wickham 2011)
rmarkdown (Allaire et al. 2018)

References

Allaire, JJ, Yihui Xie, Jonathan McPherson, Javier Luraschi, Kevin Ushey, Aron Atkins, Hadley Wickham, Joe Cheng, and Winston Chang. 2018. Rmarkdown: Dynamic Documents for r. https://CRAN.R-project.org/package=rmarkdown.

Dowle, M, A Srinivasan, T Short, S Lianoglou with contributions from R Saporta, and E Antonyan. 2015. Data.table: Extension of Data.frame. https://CRAN.R-project.org/package=data.table.

Grolemund, Garrett, and Hadley Wickham. 2011. “Dates and Times Made Easy with lubridate.” Journal of Statistical Software 40 (3): 1–25. http://www.jstatsoft.org/v40/i03/.

Müller, Kirill. 2017. Here: A Simpler Way to Find Your Files. https://CRAN.R-project.org/package=here.

———. 2018. Hms: Pretty Time of Day. https://CRAN.R-project.org/package=hms.

R Core Team. 2016. R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing. https://www.R-project.org/.

Staffell, Iain. 2017. “Measuring the Progress and Impacts of Decarbonising British Electricity.” Energy Policy 102 (March): 463–75. https://doi.org/10.1016/j.enpol.2016.12.037.

Staffell, Iain, and Stefan Pfenninger. 2018. “The Increasing Impact of Weather on Electricity Supply and Demand.” Energy 145 (February): 65–78. https://doi.org/10.1016/j.energy.2017.12.051.

Wickham, Hadley. 2009. Ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York. http://ggplot2.org.

Xie, Yihui. 2016a. Bookdown: Authoring Books and Technical Documents with R Markdown. Boca Raton, Florida: Chapman; Hall/CRC. https://github.com/rstudio/bookdown.

———. 2016b. Knitr: A General-Purpose Package for Dynamic Report Generation in r. https://CRAN.R-project.org/package=knitr.

NZ generation share trends

Electricity Authority (EA) data

Ben Anderson (dataknut@icloud.com)

Last run at: 2024-08-27 23:14:13.264296